{
 "metadata": {
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  },
  "orig_nbformat": 2,
  "kernelspec": {
   "name": "python385jvsc74a57bd04c36a90a135e4aa4256adc3692917ccc380701b533777f3da43dc8bc5264a29c",
   "display_name": "Python 3.8.5 64-bit (conda)"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2,
 "cells": [
  {
   "source": [
    "# 第二次作业第1题\n",
    "作者：赵飞\n",
    "学号：202018019427078\n",
    "日期：2021-04-24"
   ],
   "cell_type": "markdown",
   "metadata": {}
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "非线性方程组的解[a, d, r, L]依次为：\n [1.5940638047939404, 140.0, 97.86308398123771, 156.0]\n解代入后各方程误差：\n [-1.4033219031261979e-13, -7.105427357601002e-13, 0.0, 0.0]\n"
     ]
    }
   ],
   "source": [
    "from scipy.optimize import fsolve\n",
    "from math import cos, sin\n",
    "\n",
    "# 1. 直接求解该线性方程组\n",
    "def ferr(x):\n",
    "    '''\n",
    "    待求解非线性方程组的误差函数\n",
    "    '''\n",
    "    a, d, r, L = x.tolist()\n",
    "    return [\n",
    "        cos(a)-1+d**2/(2*(r**2)),\n",
    "        L-a*r,\n",
    "        d-140,\n",
    "        L-156\n",
    "    ]\n",
    "\n",
    "x1 = fsolve(ferr, [1., 1., 1., 1.])\n",
    "print(\"非线性方程组的解[a, d, r, L]依次为：\\n\", x1.tolist())\n",
    "print(\"解代入后各方程误差：\\n\", ferr(x1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "output_type": "stream",
     "name": "stdout",
     "text": [
      "(导入雅可比矩阵)非线性方程组的解为[a, d, r, L]依次为：\n [1.5940638047944053, 140.0, 97.86308398120879, 156.0]\n(导入雅可比矩阵)解代入后各方程误差：\n [-4.440892098500626e-16, -1.1368683772161603e-13, 0.0, 0.0]\n"
     ]
    }
   ],
   "source": [
    "# 2. 导入参数雅可比矩阵后再次求解\n",
    "def jacobi(x):\n",
    "    '''\n",
    "    计算非线性方程组对应的雅可比矩阵\n",
    "    '''\n",
    "    a, d, r, L = x.tolist()\n",
    "    return [\n",
    "        [-sin(a), d/(r**2), -(d**2)/(r**3), 0],\n",
    "        [-r, 0, -a, 0],\n",
    "        [0, 1, 0, 0],\n",
    "        [0, 0, 0, 1]\n",
    "    ]\n",
    "\n",
    "x2 = fsolve(ferr, [1., 1., 1., 1.], fprime=jacobi)\n",
    "print(\"(导入雅可比矩阵)非线性方程组的解为[a, d, r, L]依次为：\\n\", x2.tolist())\n",
    "print(\"(导入雅可比矩阵)解代入后各方程误差：\\n\", ferr(x2))"
   ]
  }
 ]
}